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Abstract 



?H . In this paper we consider the numerical solutions for a class of 

jump diffusions with Markovian switching. After briefly reviewing nec- 
essary notions, a new jump-adapted efficient algorithm based on the 
Euler scheme is constructed for approximating the exact solution. Un- 
der some general conditions, it is proved that the numerical solution 
through such scheme converge to the exact solution. Moreover, the 
order of the error between the numerical solution and the exact solu- 
tion is also derived. Numerical experiments are carried out to show 
^— i , the computational efficiency of the approximation. 
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> 1 Introduction 

On 

■ Stochastic hybrid systems arise in numerous applications of systems with 

£NJ ! frequent unpredictable structural changes, e.g. flexible manufacturing sys- 

| terns, air traffic management, electric power networks, mathematical finance 

and risk management etc. In typical hybrid system, the state space consists 
of two parts: one component takes values in Euclidean space while the other 
takes discrete values. For various applications of stochastic hybrid systems, 
we refer to Q31 \T7\ EB] ■ 
^ . One important classes of stochastic hybrid systems is the following jump 

diffusion with Markovian switching: 

y(t) =2/0+/ b(y(s),r(s))ds+ I a{y{s),r{s))dW{s) 
Jo Jo 

+ / / g(y(s-),r(s),v)N(ds,dv), 
Jo Jt 



(1.1) 



where r(t) is a continuous-time Markov chain represents the environment 
or modes. Such processes have been increasingly used in modelling of the 
stochastic systems which are affected by randomly occuring impulses as 
well as frequent regime switching usually modelled by Markov chain. A 
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typical application of (jl.ip stems from risk theory in insurance and finance 
where y(t) is regarded as a surplus process in Markov-modulated market. 
Applications have also been reported in a wide range of fields such as option 
pricing J3[ [19] and flexible manufacturing systems [JJ . 

Since in many problems such Levy stochastic differential equations with 
Markovian switching cannot be solved explicitly, it is important, from the- 
oretical point of view, and even more for the sake of various applications, 
to find their approximate solutions in an explicit form or in a form suit- 
able for application of some numerical methods. Although the numerical 
solutions of stochastic differential equations have been the focus of enor- 
mous research (see e.g. [Til 021 EE])) there has been significantly 
less work on the numerical methods for diffusion processes or jump diffu- 
sion processes with Markovian switching. In the literatures addressing these 
problems, Yuan and Mao |21j firstly considered the numerical solutions of 
the following stochastic differential equations with Markovian switching 

dy(t) = f(y(t),r(t))dt + g(y(t),r(t))dW(t), 
P{r(t + 5)= j\r{t) = i} = qij S + o(5) i^j. 

It was proved that the numerical solutions of (jl.2p under Euler method 
converge to the exact solutions, and the order of error was also estimated. 
Krystul and Bagchi [12] extended (|1.2p to a more general switching diffusion 
processes with state-dependent switching rates, and a approximation scheme 
for first passage time of the processes was derived. Recently, Yin, Song and 
Zhang [20] proposed a numerical algorithm for (jl.ip . and proved that the 
algorithm converges to the desired limit by means of a martingale problem 
formulation. 

In this paper, we try to develop a jump-adapted numerical algorithm for 
(jl.ip . which is based on the Euler scheme. Different from |20| . the algorithm 
we proposed does not rely on a constant-step-size scheme, instead, a series 
of jump-times driven by a Poisson random measure are taken into account. 
Another difference is that convergence properties in [2D] was proved in weak 
sense while we prove it in strong sense. Meanwhile this paper focus on the 
order of convergence, which is enlightened by the methods in [21] . 

The rest of the paper is organized as follows. Section 2 briefly reviews 
the existence and uniqueness of equation (jl.ip and describe the modified 
Euler algorithm of numerical solutions. Section 3 shows that the numerical 
solution under such algorithm converge to the exact solution in L 2 , and 
the order of the error between the numerical and exact solution is derived. 
Section 4 gives two numerical examples to examine the performance of the 
algorithm described in Section 3. 
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2 Preliminary and algorithm 



Let (Q, J-, {J 7 t}t>o, P) be a complete probability space with a filtration 
{Ft}t>o satisfying the usual conditions. Suppose that there is a finite set 
S = {1, 2, TV}, representing the possible regimes of the environment. Let 
r be a compact subset of M. not including the origin, and b(-, •) : M. d x S — > M. d , 
a(-, •) : R d x S ->• R dxd , g(; •, •) : R d x 5 x T -> M d . Consider the following 
jump diffusion with Markovian switching of the form 



!/(*) =2/o+ / 6(i/(a),r(a))(fa+ / ff(y(s),r(*))cW(a) 
Jo Jo 

+ / / p(y(a-),r(a),t;)JV(dfl,du), 

Jo Jr 



(2.1) 



with initial value y(0) = yo € R d and r(0) = io 6 5, where PF(i) = 
• • • , PF rf (i)) T is a d-dimensional J-^-adapted standard Brownian mo- 
tion, N(dt, dv) is a J-^-Poisson point process on M + x T with determinis- 
tic characteristic measure H(dv) on a compact set rcl, let N(dt, dv) = 
N(dt, dv) —U(dv)dt be the compensated Poisson random measure on R + x T, 
and r(t) is a continuous-time Markov chain taking value in a finite state 
space S with the generator Q = (qij)NxN given by 



p{ r (t + 5) =j\r(t) = i} = < 
provided 8 I 0, and 



qij5 + o(5), ifi^j, 
l + &£ + o(<5), ifi = j, 



-9m = ^Qij < +°°- 

We assume W(t),N(t, ■) and r(t) are independent. Denote iV t = N(t,T) 
by a homogeneous Poisson process with intensity A = II(r). It is well known 
that the compensated measure of N(dt, dv) is H(dv)dt, hence the centered 
Poisson measure is N(dt,dv) which is a martingale on (f2, F, {Ff}t>o, P) 
(see [TO]). 

Remarks 2.1. when d = 1, we assume that the characteristic measure 
H(dv) = XF(dv) where F{dv) is a probability distribution on V C M + and 
A is a positive constant. In this case, the integral with respect to the ran- 
dom measure N(dt, dv) is simply a compound Poisson process. We have 
J*q f r vN(dt,dv) = X^Si£i> where {N t ,t > 0} is a process with intensity 
A and G N} is a sequence of i.i.d. random variables with common 

distribution F. 
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2.1 Existence and uniqueness 

Under certain conditions we can establish the existence of a pathwise 
unique solution of (|2.ip . Here we make the following global Lipschitz and 
linear growth assumptions: 

(HI) For all (i,x,y,v) G S x M. d x R d x T, there exists a constant L\ > 
such that 

\b(x, i)-b(y, i)\ 2 + \a(x, i)-a(y, { )\ 2 +J r \d{ x i h v)-g(y, i, v)\ 2 U(dv) < L x \x-y 

(H2) For all (i, x,v) G S x R rf x F, there exists a constant L2 > such that 
J \g(x,i,v)\ 2 U(dv) < L 2 (l + |x| 2 ). 

(%3) Let ii" be the support of g(-,-,-) and U be the projection of K on 
r C R, then assume that H(U) < 00. 

We denote by | ■ | the Euclidean norm for vectors or the trace norm for 
matrices. 

Remarks 2.2. It's easy to show that if b(-, •), <r(-, •) satisfy global Lipschitz 
condition, then they also satisfy linear growth condition, i.e. for all (x, i) € 
R rf x S, there exists a constant L3 such that 

\b{x,i)\ 2 + \a{x,i)\ 2 < L 3 (l+ |x| 2 ). 

To show this, just let L 3 = 2L X V max(2|6(0, i)\ 2 + 2[cr(0,i)j 2 : i £ S). For 
convenience, we will set L = L\ V L2 V L3 in the following of this paper. 

Theorem 2.1. Ifb(x,i), a(x,i) and g(x,i,v) satisfy the conditions (HI), 
(H2), (H3), and suppose W(t),r(t), N(dt,dv) be independent. Then there 
exists a unique d- dimensional Ft- adapted right- continuous process y(t) with 
left-hand limits which satisfies equation (2.1) and such that y(0) = yo and 
r(0) = io a.s. 

Proof. See [7] □ 

2.2 Algorithm 

Now we turn our attention to numerical algorithm. Given A > as 
a step size, denote {^}«>i the usual equidistant time discretization of a 
bounded interval [0, T], i.e. t' = 0, t' i — t' i _ l = A, if t' n _ 1 < T < t' n then set 
t' n = T. Suppose r'j = inf{t > 0, N t > j},j = 1,2, •• • are the jump times 
of N t = N(t, r). Apparently, Tj = rj AT, j = 1, 2, ■■■ are also stopping 
times. Now we take a new time discretization {ifc}fc>i = {^}i>i U {Tj}j>i 
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satisfying t k+ \ — t k < A for all k > 1. This type of grid using the jump-times 
as additional grid points is usually refered to as a jump-adapted method, 
see for example [5]. 

Since N(dt,dv) and r(t) are independent, for given partition {t k } k >i, 
{r(t k ), k = 1, 2, • ■ • } is a discrete Markov chain with transition probability 
matrix (P(i,j)) NxN , here P(i,j) = P(r(t k+1 ) = j\r(t k ) = i) is the r/th 
entry of the matrix e^ tk+1 ~ tk " >( ^ , thus we could use following recursion proce- 
dure to simulate the discrete Markov chain {r(t k )} k >i, suppose r(tk) = i% 
and generate a random number £ which is uniformly distributed in [0,1], 
then we define 



r(t 



r i 2 , if i 2 € S - {N} and Y!-=i P(h,j) < £ < HhJ) 



[N, if £pP(ii, ;)<*. 

Repeating this procedure a trajectory of {r(^)}fc>i can be simulated. 
Now we could define the Euler approximation solution of (|2.ip . Let 



x to = yo, n = io, 

X tk+1 = X tk + b(X tk ,r tk )(t k+1 - t k ) + a(X tk ,r tk )(W(t k+1 ) - W(t k )) 

^ (2.2) 
+ 2^ 9 (X tk , r tk , vj )I{ tk+1=Tj } ■ 
i=i 

Let 

X(t) = X tk , r(t) = r tk (t k <t< t k+1 ). 

Hence we can define the approximation solution X(t) on the entire in- 
terval [0, T] by 



X(t) = y + f b(X(s),r(s))ds+ f a(X(s),r(s))dW(s) 
Jo Jo 

+ [ [ g(X(s-),r(s),v)N(ds,dv). 
Jo Jr 



(2.3) 



Remarks 2.3. J2f=i d(Xt k , rt k , v j)I{t k+1 =Tj} onr y take nonzero values on stop- 
ping times {tj}j>i, thus X(t) is continuous on [0, T] \ {Tj}j>\. Moreover, 
we have X(t k ) = X(t k ) = X tk . 

3 Convergence with the Lipschitz and linear growth 
conditions 

In this section, we will prove that the numerical solution X(t) converges 
to the exact solution y(t) in I? as step size A \. 0, and the order of conver- 
gence is one-half, i.e. 

E{ sup \X(t) -y(t)\ 2 ) < CA + o(A). (3.1) 

0<t<T 
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To begin with, we need the following lemma. 

Lemma 3.1. Under conditions (Til), (Ti2) and (Ti3), there exists a constant 
M which is dependent on T, L, y$, but independent of A, such that 

E( sup \y{t)\ 2 ) V E( sup \X(t)\ 2 )<M. (3.2) 

0<t<T 0<t<T 

Proof. From Holder inequality, we have 

\y(t)\ 2 = \yo+ f b(y(s),r( S ))ds + f a(y( S ),r(s))dW(s) 
Jo Jo 

+ 11 g(y(s-),r(s),v)N(ds,dv)\ 2 
Jo Jt 

<My \ 2 + 4\[ b(y(s),r(s))ds\ 2 + 4| f a(y(s),r(s))dW(s)\ 2 
Jo Jo 

+ 4| f [ g(y(s-),r(s),v)N(ds,dv)\ 2 



Jo Jr 

<4\y \ 2 +tt[ \b(y(s),r(s))\ 2 ds + 4\ [ a(y(s),r(s))dW(s)\' 
Jo Jo 

+ 4| f [ g(y(s-) 1 r(s),v)N(ds,dv)\ 2 . 
Jo Jr 

Thus for any < T' < T, we have 
\y(t)\ 2 < A\y \ 2 + AT \ 

0<t<T 



sup |y(t)| 2 <4|y | 2 + 4T / \b(y(s), r(s))\ 2 ds 
<t<T' Jo 

+ 4 sup | / cr(y(s),r(s))dW(s)\' 

0<t<T' JO 



(3.3) 



o<t<T' Jo Jr 
From (Til) and Remark 2.2, we have 

rT' 
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ft 

+ 4 sup | / / g(y(s-),r(s),v)N(ds,dv)\ 2 



E [ \b(y(s),r(s))\ 2 ds <E f L(l + \y(s)\ 2 )ds 
Jo Jo 

< Ci + Ci / E\y(s)\ 2 ds. 
Jo 



(3.4) 



Here C\ is a positive constant, in fact, C\ = T'LVL in (|3.4p . For convenience, 
in the following of this paper, we may frequently denote the related constants 
by Cfc. In that case, we only mean there exist such a positive constant rather 
than some specific one. 
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Since / a(y(s), r(s))dW(s) is martingale (see [6]), therefore, from Doob 
Jo 

martingale inequality and ("HI), we have 



E{ sup | / <r(y(s),r(s))dW(s)\ 2 ) 

0<t<T' JO 

<AE\ / a(y(s),r( S ))dW(s)\ 2 
Jo 

= AE / \a(y(s),r(s))\ 2 ds (3.5) 
Jo 



<AE ! L{l + \y(s)\ 2 )ds 
Jo 

<C 2 + C 2 E\y(s)\ 2 ds. 
Jo 



For the last term of (|3,3|) . notice that N(dt, dv) = N(dt, dv) — U(dv)dt is 

the compensated Poisson random measure, and / / g(y(s~),r(s),v)N(ds,dv) 

Jo Jr 

is martingale (see we have 



E{ sup | / / g(y( 
o<t<T' Jo Jr 



S -),r(s),v)N(ds,dv)\ 2 ) 



= E( sup | / / g(y(s ),r(s),v)N(ds,dv) 
o<t<T' Jo Jr 

+ 11 g(y(s-),r(s),v)n(dv)ds\ 2 ) (3.6) 
Jo Jr 

<2E{ sup | [ [ g(y(s~),r(s),v)N(ds,dv)\ 2 ) 
o<t<T' Jo Jr 



+ 2E( sup | / /■ 

o<t<T> Jo Jr 



g(y(s ),r(s),v)U(dv)ds\ ). 

Hence, by Doob martingale inequality, we have 

E{ sup | / / g(y(s-),r(s),v)N(ds,dv)\ 2 ) 
o<t<T' Jo Jr 

rj-if 

<4E\ [ [ g(y( S -),r(s),v)N(ds,dv)\ 2 (3.7) 
io Jr 

rpl 

= AE J J \g(y(s-),r(s),v)\ 2 U(dv)ds. 
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On the other hand, we can get 

E{ sup | f [ g(y(s-),r(s),v)n(dv)ds\ 2 ) 

0<t<T' Jo Jr 

<E( sup tf f \g(y(s-),r( S ),v)\ 2 U(dv)ds) (3.8) 
o<t<T' Jo Jr 

<TE [ [ \g(y(s-),r(s),v)\ 2 n(dv)ds. 
Jo Jr 

Combining fl3jl, ([321), $M and C^ 2 ), we obtain that 



E( sup | f [ g{y(s-),r(s),v)N{ds,dv)\< 
o<t<T' Jo Jr 



T' 

2i 



<(8 + 2T)E / \g(y(s-),r(s),v)\ 2 n(dv)ds 
Jo Jr 



(3.9) 



< (8 + 2T)E [ L(l + \y(s-)\) 2 ds 
Jo 

<C 3 + C 3 E\y(s-)\ 2 ds. 
Jo 

Take expectation of (|3.3p and substitute (J33]), (|3.5|) . (|3.9|) into it, we have 



E( sup \y(t)\ 2 ) 

0<t<T' 



< 4|y | 2 + (Ci + C 2 ) J E\y(t)\ 2 dt + C 3 J E\y(t-)\ 2 dt (3 . 10) 
<Ci + C 4 E{ sup |y(s)| 2 )cft. 

Jo 0<s<t 

From Gronwall inequality, we can show that 

E{ sup \y{t)\ 2 ) < C 4 (l + OTe C4T ) = M x . 

0<t<T 

Similar to above procedure, we know that for X(t) 

E{ sup \X{t)\ 2 ) <4|y | 2 + (d + C 2 ) F E\X(t)\ 2 dt + C 3 C E\X{t~)\ 2 dt. 

0<t<T' Jo Jo 

According to Remark 2.3 

X(t) = X(t~) =X tk =X(t k ), t k <t<t k+1 , 

X(t) = x tk = x{t k ), x(t-) = x tk _ x = X(t k -i), t = t k . 



In both cases we have 

\X(t)\ 2 V \X(t-)\ 2 < sup \X(s)\ 2 . 

0<s<t 

Thus by Gronwall inequality 

E( sup \X(t)\ 2 ) < M 2 . 

0<t<T 

Finally there exists a positive constant M such that 

E{ sup \y(t)\ 2 )VE( sup \X(t)\ 2 )<M. 

0<t<T 0<t<T 

□ 

Now we come back to our main theorem, the proof below mainly refers 
to Yuan and Mao [21j . 

Theorem 3.1. Assume the jump diffusion with Markovian switching \2. 1\) 
defined on {£l,F,{Ft}t>o,P) satisfying 

a) W(t),r(t), N(dt,dv) are independent, 

b) &(•,•), g(-,-,-) satisfy conditions (HI), (T~L2) and (H3), 

then the unique strong solution y(t) and numerical solution X(t) obtained 
in section 2.2 satisfying: 

E( sup \X{t) -y(t)\ 2 ) < CA + o(A). 

0<t<T 

where C is a positive constant independent of A. 

Proof. For every < T' < T, similar to Lemma 3.1 we can easily verify that 
E{ sup \X{t) - y{t)\ 2 ) 

0<t<T' 

< 3TE / \b(X(s),r(s)) - b(y(s),r(s))\ 2 ds 

r (3.11) 

+ 12E / \a(X( S ),r(s))-a(y(s),r( S ))\ 2 ds 
Jo 

+ CE [ [ \g(X(s-),r(s),v)-g(y(s~),r(s),v)\ 2 U(dv)d S . 

Jo Jt 
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We focus on the third part of the right side 



E [ [ \g(X(s~),f( S ),v)-g(y(s~),r(s),v)\ 2 U(dv)d S 
Jo Jr 

< 2E J o J r \g(X(s-),r(s),v) - g(y( S -),r( S ),v)\ 2 U(dv)ds (3.12) 
+ 2E J o J r \g(X(s-),f(s),v) - g(X(s-),r(s),v)\ 2 U(dv)ds. 
By condition (Til), we have 

rp! 

E [ [ \g(X(s-),r(s),v) - g(y( S ~),r(s),v)\ 2 U(dv)d S 
Jo Jr 

rT' 



< 2L 2 f E\X(s~) -y(s-)\ 2 ds. 
Jo 



(3.13) 



By (U2), we can compute that 

r-T 



zJ Q J r \g(X(s-),f(s),v) - g(X(s'),r(s),v)\ 2 U(dv)ds 

= J\g(X( 8 -),f( 8 ),v) 

-g(X(s-),r(s),v)\ 2 U(dv)ds | N T , n, ...t Nt }] 
= E[Y^E[f k+1 f \g(X(s-),f(s),v) 



fc>0 Jtk 

|2i 



-g(X(s ),r(s),v)\ 2 U(dv)ds \ N T , n, ...t Nt }] 

^ /■**+! r (3-14) 
= E[J2E[ \g(X(t k ),r(t k ),v)-g(X(t k ),r( S ),v)\ 2 

k>0 ^ tk ^ r 

x /{ r ( s )^ r (i fc ) } n((iw)ds | N T ,n,...T NT ]] 

< 2E\Y^E[ f k+1 f (\g(X(t k ),r(t k ),v)\ 2 + \g(X(t k ),r(s),v)\ 2 ) 

k>0 tk r 

x I {r{s )^ r{tk )}Ii{dv)ds | N t ,ti,...t Nt \] 

< CE[J2E[ / [1 + \X(t k )\ 2 }I {r(s ^ r{tk)} ds I AT T ,n, ...tat t ]]. 

fc>0 tk 

Note that given {Nt, t\, ...t/v t }, 1 + |X(t^)| 2 and I^ r ^ s ^ r ^ tk ^ are condition- 
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ally independent with respect to r(t k ), thus 



E [ h+ \l + \X(t k )\ 2 ]I {r{s)j , r{tk)} ds 

Jtu 



■fc+1 



' E[E[(1 + \X(t k )\ 2 )I {r{s)7 , r(tk)} | r(t k )]]ds (3.15) 

tk+1 



tk 



£[£[(1 + |X(t fc )r) | r(t fc )]^[/{rW/r(t k )} I r (*fe)]]d 



s. 



On the other hand, by the Markov property, we get 

E[I {r{s) ^ r{tk)} | r(t k )]] =J2 I {r(t k )=i} P ( r ( s ) ¥= i I T ^k) = i) 

ies (3.16) 
< CA + o(A). 



Substituting (|3.16p into (|3.15|) gives 

"tk+i 



ftk+l 

E [l + \X(t k )\ 2 ]I {r{s)Mtk)} ds 

< [ ^ E[E[(1 + \X{t k )\ 2 ) | r(t k )](CA + o(A))]ds 

Jtu 



ftk+l 

< (CA + o(A)) / E[E[(1 + \X(t k )\ 2 ) | r(t fc )]]da 

ftk + l 

= (CA + o(A)) / S[l + |X(t fc )| 2 ]ds. 

By Lemma 3.1 there exists a positive constant M such that 
E[l + |X(t fe )| 2 ] < 1 + E( sup \X(t)\ 2 ) < M, 

0<t<T 

then 



£ / [1 + |X(t fc )| 2 ]I{r(,)^(t fc )}^ < - *jfe)(CA + o(A)). (3.17) 
't k 
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Substitute flSTZP into (1X1^) 



E / / \g(X(s-),f(s),v) - g(X(s-),r(s),v)\'U(dv)ds 

JO JT 

< CE[Y^ E[(t k +i - t k )(CA + o(A)) | N T , n , ...t Nt }\ 

k>o (3.18) 

= CE$2(t k+1 -t k )(CA + o(A))] 

k>0 

< CA + o(A). 

Here we use the fact that (tk+i — tk)(C A+o(A)) are measurable with respect 

to the o"-algebra generated by {Nt,ti, ...tn t }- 

Now substitute (pTTH]) . (I3~T8]1 into (ETT21 . we obtain that 



/ |sr(X(s-),r( S ),^) - 5 r( 2/ ( S -),r(s),t;)| 2 n( C ?-u)d S 

jo 

< 2L 2 / £|X(0 - y{s~)\ 2 ds + CA + o(A). 

Similar to (|3.19p . we have 
-T' 
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E I \b(X(s),r(s))-b(y(s),r(s))\ 2 ds 



< 2L 2 [ E\X(s) - y{s)\ 2 ds + CA + o(A), 
Jo 

E / \a(X( S ),r(s))-a(y(s),r(s))\ 2 ds 



o 



< 2L 2 [ E\X(s) - y(s)\ 2 ds + CA + o(A). 
Jo 

Substituting ([319]) . (13301) . (pT211 into (I3TTTD shows that 



(3.19) 



(3.20) 



(3.21) 



E( sup \X(t)-y(t)\ 2 ) 

0<t<T' 

r T ' 

<C E\X(s)-y(s)\ 2 ds (3.22) 
Jo 

+ C / £|X(s~) - y(s-)\ 2 ds + CA + o(A). 

Note that 

E\X(s) - y(s)\ 2 < 2E\X(s) - y(s)\ 2 + 2E\X(s) - X(s)\ 2 , (3.23) 
E\X(s-) - y(s-)\ 2 < 2E\X(s-) - y( S -)\ 2 + 2E\X{ 8 ~) - X( S -)\ 2 . (3.24) 

12 



Suppose t k < s < tk+i, then 
\X(s)-X(s)\ 2 = \X(s)-X(t k )\ 2 



b(X{u),f{u))du+ / a(X(u),f(u))dW(u) 

+ f f g(X(u~),f(u),v)N(du,dv)\ 2 
Jt k Jr 

= \KX tk ,r tk )(8 - t k ) + a(X tk ,r tk )(W(s) - W(t k )) 

+ [ [ g{X tk ,r tk ,v)N{du,dv)\ 2 
Jt k Jr 

< 3\b(X tk ,r tk )\ 2 (s - t k ) 2 + 3\a(X tk ,r tk )\ 2 (W(s) - W(t k )f 

+ 3| f / g(X tk ,r tk ,v)N(du,dv)\ 2 . 
Jt k Jr 



Since for given {t k }k>i, W(s) — W(t k ) and (X(t k ),rt k ) are independent, we 
have 

E\X(s)-X(s)\ 2 

< 3E[\b(X tk ,r tk )\ 2 (s - t k ) 2 + \a(X tk ,r tk )\ 2 (W( S ) - W(t k )) 2 



+ I f [ g(X tk ,r tk ,v)N(du,dv)\ 2 ] 
Jt k Jr 

3^[|6(X tfc ,r t J| 2 ( S -i fe ) 2 + KX tfc ,r t J| 2 (W( S )-W(^ 



+ | / / 0(X tfc ,r tfc ,w)JV(dMw)| 2 |i fe ]] 
Jt k Jr 

< 3E[E[\b(X tk ,r tk )\ 2 A 2 \ t k ] + E[\a(X tk ,r tk )\ 2 (W(s) - W(t k )) 2 \ t k \ 

+ E[\ f [ g(X tk ,r tk ,v)N(du,dv)\ 2 \t k ]] 
Jt k Jr 

= 3A 2 E[\b(X tk ,r tk )\ 2 ] + 3E[EMX tk ,r tk )\ 2 | t k ]E[(W(s) - W{t k )) 2 \ t k ]} 

+ 3E[E[[ S [ \g(X tk ,r tk ,v)\ 2 U(dv)du\t k ]} 
Jt k Jr 

= 3A 2 E[\b(X tk ,r tk )\ 2 ] + 3E[E[\a(X tk ,r tk )\ 2 \ t k ](s - t k )} 
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+ 3E[E[J^\g(X tk ,r tk ,v)\ 2 U(dv) \ t k ](s - t k )] 
<3E[\b(X tk ,r tk )\ 2 ] + 3E[E[\a(X tk ,r tk )\ 2 \t k ]A] 

+ 3E[E[^\g(X tk ,r tk ,v)\ 2 U(dv) | t k ]A] 
= 3A 2 E[\b(X tk ,r tk )\ 2 ] + 3AE[\a(X tk ,r tk )\ 2 } 

+ 3AE[J^\g(X tk ,r tkJ v)\ 2 U(dv)}. 

Hence by (Til) and Lemma 3.1, we get 

E\X(s) -X(s)\ 2 < CA + o{A). (3.25) 
Similarly, we can show that 

E\X(s~) -X(s~)\ 2 < CA + o(A). (3.26) 

Substituting §£23$ , §£M$ , (|5^5]> . (13T2B1 into (l3T22"j) immediately shows 
that 

E( sup \X{t) - y{t)\ 2 ) 

0<t<T' 

r T ' 

< C / (E\X(s) - y(s)\ 2 + E\X(s~) - y{ S -)\ 2 )ds + CA + o(A) 
Jo 

<C E{ sup \X(s) - y(s)\ 2 )dt + CA + o(A). 

Jo 0<s<t 

Therefore, from Gronwall inequality we obtain that 

E( sup \X(t) -y(t)\ 2 ) < CA + o(A). 

0<t<T 

The proof is complete. □ 



4 Numerical examples 

In this section, we discuss two numerical examples to illustrate our the- 
ory established in the previous sections. It is shown that the computer 
simulation based on our numerical method is feasible and efficient. 



14 



Example 4.1. Consider the following one-dimensional geometric Levy pro- 
cesses 



y(t) = Vo+ y(s)fi(r(s))ds + / y(s)a(r{s))dW{s) 

Jo Jo (41) 

+ / y(s-)g(r(s))dN(s), 
Jo 

with initial value y(0) = yo, r(0) = io, where W{t) is a scalar Brownian 
motion, N(t) is a Poisson process with intensity A > and r(t) is a right- 
continuous Markov chain taking values in S, suppose that r(t), W(t), N(t) 
are independent. It is obvious that the assumptions ("HI), (%2), ("H3) are 
satisfied for equation (|4.ip . It is well known that equation (|4.ip has an 
explicit solution 

1 ' 2/ 



y(t)=y exp[ fi(r{s))ds + a(r(s))dW(s) - - a 2 (r(s))ds] 
Jo Jo z Jo 

N{t) (4.2) 



n(l + 5 (r(x„))), 



X 

n=l 



where r n are jump times determined by Poisson process N(t). 

As described in Section 2, we get the numerical solution of (|4.ip 



= 2/0 + y X (a)M(f (s))ds + / X(s)a(f(s))dW(s) 
X(s)g(f(s))dN(s). 



(4.3) 



By virtue of Theorem 3.1, it follows that the solutions y(t) and X(t) are 

close in the L 2 -norm, with the rate less than CA when A J, 0. To examine 

this convergence result, we firstly simulate the trajectory of y(t) and X{t) for 

sufficient small stepsize A > 0, then use simulated value E{ sup \X(tk) — 

t k e[o,T] 

y{tk)\ 2 ) a t new division points tj, to estimate E( sup \X(t) — y(t)\ 2 ), where 

0<t<T 

E denotes sample mean. 

For simulation reason, it is convenient to transform (|4.2|) into following 
recursion form with y(to) = yo, 



Z/0fc+i) ~ y(t k )exp[(t k+1 - t k )n(r tk ) + (W(tfc+i) - W{t k ))a(r tk ) 

N(T) 

- ~(t k+1 - t k )a 2 {r tk )} x (1 + 9{n k+1 ))I{t k+1 = Tj). 



Notice that y(tfc + i) in (|4.4p is not the exact value of y(t) at the division 
points t k+ i, because r(s) is not necessarily constant on [t k ,t k+ i\. However, 
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since 



P{r(t k+1 ) = i\r(t k ) =i} = l + qu(t k+1 - t k ) + o(t k+1 - t k ) -)• 1, 

as A \, 0, for sufficiently small A, it is reasonable to use (|4.4p as an approx- 
imation of the exact solution of y(t). 

For the discrete approximate solution (|4.3p at the division points t k+ \, 
we have 



X tk+1 = X tk + X tkl x{r tk ){t k+1 - t k ) + X tk a(r tk )(W(t k+1 ) - W(t k )) 

N(T) 

+ S X tk9(rt k )I{t k+1=Tj }- 

3=1 



(4.5) 



Now we take specific data to computer the above example. Let r(t) be 
a right-continuous Markov chain taking values in 5 = {1,2} with generator 



Q 



-0.5 0.5 
0.5 -0.5 



therefore, the one step transition probability from r(t k ) to r(t k+ i) is e^ tk+1 * & ) < 3. 
Let /j,(r(t)), a(r(t)), g{r{t)) take values as follows 



//(l) = 0.15, a(l) = 0.1, <?(!) 



-0.2; 



/i(2) = 0.05, <j(2)=0.1, 5(2) = -0.1. 

Let the Poisson intensity A = 1, and choose initial values yo = 10, 
ro = 1, T is fixed at 10. By applying the previously described procedure, in 
Figure 1 the trajectory of the approximate solution X(t) with A = 0.01 is 
constructed. 




Figure 1. The trajectory of the approximate solution when A = 0.01 
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To carry out the numerical simulation we successively choose the stepsize 
A as the following Table 1, and for each A, we repeatedly simulate and 

compute sup (\X(tk) — y{tk)\ 2 ) for 1000 times, then calculate the sample 

t fc e[o,io] 

mean E( sup \X(tk) — y(tk)\ 2 ). The results are listed in the following 

* fc e[o,io] 

Table 1. 



A 


E( sup \X{t k ) - y{t k )\ 2 ) 

t fc 6[0,10] 


0.001 


0.000179654 


0.005 


0.000908019 


0.01 


0.001793076 


0.02 


0.003626408 


0.03 


0.005427134 


0.05 


0.009155912 


0.08 


0.014817250 


0.1 


0.018204170 



Table 1. Estimation of the error between the numerical and exact solutions 



Clearly, the numerical method reveals that the numerical solution X(t) 
converges to the exact solution y{t) in I? as step size A \. 0, and the order 
of convergence is one-half, i.e., 

E{ sup \X(t) -y{t)\ 2 ) < CA + o(A), 

0<t<T 

which strongly demonstrate our result. 

Example 4.2. In this example, as an application, we turn our attention to 
the expected ruin time for a surplus process with regime switching. Denote 
by r(t) the external environment process, which influences the frequency of 
claims and the distribution of claims. Suppose that the process r(t) is a ho- 
mogeneous, irreducible and recurrent Markov process with finite state space 
S = {1, 2, ■ ■ ■ , iV} with generator Q = {qij) NxN and stationary distribution 
7r = (7Ti,7T2,--- , 7Tjv). N(t) is a Markov- modulated Poisson process which 
accounts the arrival of claims, i.e. at time t, claim arrival occur accord to 
a Poisson process with constant intensity Aj when r(t) = i, and the corre- 
sponding claim size Uj have distribution -Fj(x) with means /Xj. The surplus 
process is given by 

N(t) 

Y(t) = u + t-^2u n , (4.6) 

n=l 

where the aggregate premium received with rate 1 during interval (0, t] and 
u is the initial reserve, the i.i.d. random variables U n and N(t) are condi- 
tionally independent given r(t). 
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Define T* = inf{t > | Y(i) < 0} to be the time of ruin of (gSJ), and 
define the expected ruin time, given that the initial environment state is 
i and the initial reserve is u, by £,i(u) = E[T* | r(0) = i,Y(0) = u]. Let 
i] = YliLi ^iXi^i and p = f]^ 1 — 1. It's known that p < ensures Ci( u ) < °°- 

Now we consider the special case of two environmental states. Let 

Q=[ - qi qi 

then the stationary distribution of r(t) is tt\ = — — — , 1T2 = — — — Sup- 

qi + q2 qi+ q2 

pose that claim sizes U\,U2,-" are independent and exponentially dis- 
tributed with mean p, which are independent of the environment process 
for convenience. Thus according to [3], we have 

Proposition 4.1. If p < 0, then it holds 

V ~ (47) 
Uu) = A 2 + ^- + BD(k)e hu , 

Tj — 1 

where k is the unique negative root of the following equation 

P(k) = k 3 + k 2 {pi +p 2 -q\- q 2 ) 

+ k{pip 2 - P\q2 ~ P2q\ ) (4.8) 

p p p p 

= 0, 

with pi = Aj, i = 1,2, and (Ax,B,A 2 ) is the unique solution of the 

It- 
following liner equation 



1 



'qi 

1 1/(^ + 1) J I B I =-^ T I p I , (4.9) 
D(k)/(kp + l) 





where 



D(k) = gi + ^i+^i-l)-fcV 
qi+kpqi 



From Proposition 4.1, it is easily to calculate £i(it) and £ 2 {u) if the 
appropriate values of parameters in the model are given. However, in this 
paper we will use the numerical scheme described in Section 2 to generate 
the simulating paths and then derive the simulated value of £1 (u) and £2 (u) ■ 
Generally speaking, this simulation approach is feasible and efficient for 
obtaining numerical solutions to problems which are too complicated to 
solve analytically. Here we take for calculating £i(u) as an example. 
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Let 




(4.11) 



and Ai = 1, A2 = 2, /i = 1, then it is easy to obtain that it = (-, -), 77 = -, 

p = — — < 0, pi = 0, P2 = — 1. Solving equations (|4.8|) and (|4.9|) respectively, 
we get k = -0.6751309, A x = 2.712742, B = 1.481194. Substituting into 
(4.7) gives the exact value of £i(u) 

= 2.712742 + 2u + 1.481194 e - a6751309u . (4.12) 

For comparing the exact value (|4.12p with the simulated value of £i(it), 
we adopt the numerical algorithm in Section 2. After simulating the tra- 
jectory of Y(t), which is also denoted by X(t), it is easy to explore the 
trajectory behavior on [0, T]. To avoid unnecessary error, we choose a suf- 
ficient large interval [0, T] and record the ruin time T* of X(t) for every 
trajectory. Set T* = T if T* > T. In the following, let the stepsize A = 0.01 
and T = 100, for each initial reserve u = 5,8, 10, 15, 20, we repeatedly simu- 
late and compute T* for 1000 times respectively, then calculate the sample 
mean £i(u). The results are listed in the following Table 2. 



u 




&(«) 


5 


12.76339 


12.4102 


8 


18.71942 


18.7526 


10 


22.71447 


23.1208 


15 


32.71280 


32.5327 


20 


42.71274 


41.9659 



Table 2. Exact and simulated value of the expected ruin time for initial reserve u 



Clearly, the simulated value £i(u) is considerably close to £i(u), which 
again demonstrate the efficiency of our numerical method. 
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